A complete NLO calculation of the J j'^ and production at hadron colliders* 
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A complete next-to-leading order (NLO) calculation in a, for the J/^p and ^p' prompt production 
at the Tevatron, LHC, and RHIC in nonrelativistic QCD is presented. We argue that the next- 
to-next-to-leading order (NNLO) color-singlet contribution may not be so important as to resolve 
the large discrepancy between theory and experiment in J/ip large pr production cross sections. 
Therefore, a complete NLO calculation, including both color-singlet and color-octet contribution, is 
necessary and essential to give a good description for J/i/) and t/i' production. We also study the 
methods to fit the long-distance matrix elements using either two linear combined matrix elements 
or three matrix elements, and find these two methods can give consistent results. Compared with the 
measurements at the LHC and RHIC for prompt J/tp and tp' production, our predictions are found 
to agree with all data. In particular, the recently released large pr data (up to 60-70 GeV) at the 
LHC are in good agreement with our predictions. Our results imply that the universality of color- 
octet matrix elements holds approximately in charmonium hadroproduction, when one uses fixed 
order perturbative calculation to describe data (the data in small pt region are not included). Our 
work may provide a new test for the universality of color-octet matrix elements, and the color-octet 
mechanism in general. 

PACS numbers: 12.38. Bx, 13.85.Ni, 14.40.Pq 



I. INTRODUCTION 

Heavy quarkonium is a multiscale system which can 
probe various regimes of QCD. Thus, an understanding 
of heavy quarkonium production is particularly interest- 
ing. To solve the large discrepancy between CDF data at 
the Fermilab TevatronQ of tp' production at high pt and 
theoretical predictions, the color-octet (CO) mechanism 
3 was proposed based on nonrelativistic QCD (NRQCD) 
factorization With the CO mechanism, QQ pairs can 



be produced at short distances in CO (Sq, S[, Pj) 
states and subsequently evolve into physical quarkonia 
by nonperturbative emission of soft gluons. It can be 
verified that the partonic differential cross sections at 
leading-order (LO) in as behave as for '^S^^\ and 

1/Px for ^Sj^^ and ^Pf\ all of which decrease at high px 
much slower than 1/p^ of the color-singlet (CS) state. 
The CO mechanism could give a natural explanation for 
the observed pt distributions and large production rates 
of tp' and J/ip [sj, l6|. However, CO mechanism seems 
to encounter difficulties when the polarization of J/4> is 
also taken into consideration 0, [8|. To exploit the un- 
derlying physics, lots of efforts have been made, either 
by introducing new channels P-fTlj or by proposing other 
mechanisms [13, UM ■ 



*An expanded version based on Ref. [J]. 



It is a significant step to work out the next-to-leading 
order (NLO) QCD correction for the CS channel, which 
enhances the differential cross section by about 2 orders 
of magnitude at high PT[i3|, and changes the J/tp polar- 
ization from being transverse at LO into longitudinal at 
NL0[15]. Although the CS NLO cross section still lies 
far below the experimental data, it implies that, com- 
pared to the as suppression, kinematic enhancement at 
high Pt is more important in the current issue. This ob- 
servation is also supported by our recent workflG^ for Xc 
production, where we find the ratio of production rates of 
Xci I Xa\ can be dramatically altered by the NLO contri- 
bution due to change of the pT distribution from 1 lp\ at 
LO to 1/p^ at NLO in the CS P-wave channels. So we 
can conclude nothing definite until all important chan- 
nels in 1/pt expansion are presented. It means the CO 
channels ^Sg*' [l^ and should be considered at NLO, 

while the CS channel "^S^^ at next-to- next-to-leading or- 
der (NNLO) in as- Among these corrections, the com- 
plete NNLO calculation for the CS channel is currently 
beyond the state of the art, and instead, the NNLO* 
method is proposed fT9!,l20j. Compared to NLO, the only 
potentially not suppressed contribution within NNLO CS 
channel is gluon fragmentation, which gives a new scaling 
behavior of \lp\. for the cross section. But, as studied in 
ref. 21,], these fragmentation contributions are ignorable, 
compared with experimental prompt production data of 
J/V', and we will further argue about this point in Sec. 
IIII Al As a result, we expect a complete NLO calcula- 
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TABLE I: Experimental conditions with various experimental 
collaborations, ro and ri are theoretical predictions related 
to the short-distance coefficients. 



tion of J/ ■0 production is necessary and sufBcient to give 
a reasonable description of the experiment data. 

Currently, while J/tJj production in two-photon colli- 
sions at CERN LEP2[22] and photoproduction at DESY 
HERA [23h25i] are shown to favor the presence of CO 
contribution, the J/V' production at B factories is de- 
scribed well using NLO CS model and leaves little room 
for the CO contributions |26j-i29j. J/tp production in as- 
sociation with a VF-boson or Z°-boson at the LHC is 
also studied [30l|. However, in all previous works for 
heavy quarkonium production, CO long-distance matrix 
elements (LDMEs) were extracted at LO, which surfer 
from large uncertainties. In order to further test the CO 
mechanism, it is necessary to extract CO LDMEs at NLO 
level. This was studied in our recent work Ref. for Xcj 
and Refs.(ll.il7i] for J/-0 and ^p' . Based on Ref. Jj], we fur- 
ther study J/ ■0 and -0' hadron production including more 
detailed discussions in this work. 

The remainder of this paper is organized as follows. 
In Sec. im we perform a fit to the CO LDMEs for 
and JM using the pr distributions measured by CDF in 
Ref.[3l| and Ref. 32 1 respectively. In the fit of J/ 4', feed- 
down contributions from Xcj and '0' are considered. We 
refer interested readers to Ref.[l3l for details on the cal- 
culation and the input parameters. We will study further 
theoretical uncertainties in Sec. IIIII Then, we compare 
our predictions with new LHC data and RHIC data in 
Sec|Vl After that, a related work of NLO correction to 
J/ip production is compared with ours. We finally give 
a brief summary in Sec. IVII 



II. FIT COLOR-OCTET MATRIX ELEMENTS 
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TABLE II: Fitted Color-Octet LDMEs in J/i>{tp') produc- 
tion with p^^ = 7 GeV. Here ro = 3.9, ri = --0.56 are 
determined from short-distance coefficient decomposition at 
Tevatron. The first errors are due to renormalization and 
factorization scale dependence, while the second errors come 
from the fit. Color-Singlet (^sF) LDMEs {O") are estimated 
using a potential model result pj3l| . 



longer the unique source for the high px contribution. In 
fact, the following decomposition for the short-distance 
coefficients holds within an error of a few percent: 



(1) 



where we find tq — 3.9 and ri — —0.56 for the exper- 
imental condition with CDF at the Tevatron. ro.i for 
other conditions discussed in this work can be found in 
Table |T1 As a result, it is convenient to use two linearly 
combined LDMEs 



mi 



(2) 



when comparing theoretical predictions with experimen- 
tal data for production rates at the Tevatron and LHC. 
As pointed out in Ref. [13, ahhough both (©"'/'/'(^'f')) 
and d(T[ Pj ] depend on the renormalization scheme and 

the renormalization scale 
dent of them. 



m(^^ is almost indepen- 



We note that the curvature of experimental cross sec- 
tion is positive at large pT but negative at small pT, with 
a turning point at pT ~ 6 GeV. But the theoretical cur- 
vature is always positive. This implies that data below 
7 GeV can not be well explained in fixed order pertur- 
bative QCD calculations. If including these data in the 
fit, it will cause a large x^: which indicates the fit is not 
reliable. Therefore, in our fit we introduce a p™' and 
only use experimental data for the region pT > p^"'. In 
the following we use p^*' = 7 GeV. 

By fitting the px distributions of promptj/'' and J/ 
production measured at the Tevatron[3l|, 132] in Fig. [T] 
and Fig. [21 the CO LDMEs are determined as showing in 
Table HIl while the CS LDMEs are estimated using a po- 
tential model result of the wavefunctions at the origin[33|. 
In Fig. [Hand Fig.[2]we also give the predictions of prompt 
-0' and J/'0 production at LHC with = 14 TeV and 

|y| <3. 



3 



Tevatron 
M] Tevatron 
Mo Tevatron 
Total Tevatron 
Total LHC 
CDF Data 



3. 

T 




FIG. 1: (Color online.) Transverse momentum distributions 
of prompt tp' production at the Tevatron and LHC. CDF data 
are taken from Ref. ^31]. The LHC prediction corresponds to 
^ = 14 TeV and \yj/^\ < 3. 
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FIG. 2: (Color online.) Transverse momentum distributions 
of prompt J/t/) production at the Tevatron and LHC. CDF 
data are taken from Ref . [s^ . The LHC prediction corre- 
sponds to \/S = 14 TeV and \y^>\ < 3. 



III. THEORETIC UNCERTAINTIES 

A. Uncertainty from NNLO color-singlet 
contribution 

Ordinarily, errors come from higher order contribu- 
tions can be estimated by varying renormalization scale 
and factorization scale. This is the case for CO contribu- 



tions which have been considered in the fit. However, for 
CS contribution, new kinematic enhanced channels will 
open at NNLO which behavior as l/p^- Because the new 
channels have different px behavior from LO and NLO 
contributions, its influence can not be simply estimated 
just by varying parameters at NLO calculation. 

A complete NNLO calculation for CS is currently far 
beyond the state of the art, instead, a NNLO"^ method 
is proposed[l^[23|, in which only tree level diagrams are 
considered and an infrared cutoff (s™'") is imposed to 
control soft and coUinear divergences. As behav- 
ior channels are presented for the first time at NNLO, 
their contributions do not have divergences and should 
almost not dependent on s™™ supposing s™™ is sufh- 
ciently small. Generally, for small s™'" and large px, the 
NNLO* contributions can be expanded as 



1 V(p|/s--) 

d^NNLO* — + Cq p 

Pt Pt 



(3) 



where ... represents remained contributions which are 
not important. To demonstrate terms other than l/p^ 
have negligible contributions, authors in Ref. (l9| vary 
the s™™ and show that the yield dtTNNLO* becomes in- 
sensitive to the value of s™™ as px increases. The 
NNLO* contributions are then concluded to be large and 
important [li,!!^]. 

In the following, however, we will argue that the NNLO 
CS contribution should not be so large as the NNLO* 
method expected. We first point out that, there could be 
a misunderstanding in Ref. [19] when trying to demon- 
strate the 1/px term is the most important one. In fact, 
even if the second term in Eq. ([3]) is much larger then the 
first term, dcNNLO 
at large px, the reason is 



will also become insensitive to s"J'" 
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Thus it is needed to restudy which term is dominant in 
dcTNNLO* in the current experimental px region. 

Our strategy to study this problem is fitting the px 
behavior of 



R* = dcrNNLO*/dcrNLO, 



(5) 



where daNLO is well known to behave as 1/py at large 
px- If C4 term is dominant in daNNLO*, R* will behave 
while if cg term is dominant, R* will behave as 



2 

as Px ! 

log^(p|,/s™'"). As there is no difference between J /^p and 



T except a mass scale change, we will use the daNNLO* 
results for T in Ref. 19]. Specifically, we define 
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FIG. 3: (Color online.) Transverse momentum distributions 
of functions of /i, /2 and /a. It implies the NNLO* result is 
dominated by the double logarithm enhancement, which will 
be canceled in a complete NNLO calculation. See text for 
definition of fi. 



while /i and /2 correspond to /i and J2 normalized at 
Pt = 30 GeV. The transverse momentum distributions 
of /i and /2 are presented in Fig. [31 where we find /2 
is almost fixed to 1 when pT > 15 GeV while /i still 
varies significantly in this pt region. As a result, the 
R* behaves similar to log^(p^/s™'"). To further test this 
double logarithm behavior, we define 



h 



R^/\0g\pl/sf^)l...^2ml 
i?Vlog'(p|/s-'")L»i„=0.5m? 



(7) 



It can be found in Fig. [3] that /s is very close to when 
Pt > 12 GeV, which confirms our expectation for the 
double logarithm behavior. 

Based on the above discussion, we may conclude that 
dcTNNLO* in the current experimental pt region is dom- 
inated by cg term which has double logarithm enhance- 
ment relative to NLO result ^. The double logarithm, 
originated from IR cutoff, will be canceled in a complete 
NNLO calculation with both real and virtual corrections 
taken into consideration. Therefore, a complete NNLO 
result should have no large enhancement relative to NLO 
result 17], considering the suppression due to an extra as 
in NNLO. In other words, the NNLO* method may have 
overestimated the NNLO contributions. 

Having found that the NNLO CS contribution should 
not be large relative to the NLO contribution, we may 
ignore the theoretical uncertainty from NNLO because 
the CS NLO result is smaller than experimental data by 
at least a factor of 10 at > 7 GeV. 



B. Uncertainty from decomposing P-wave channels 

There are two reasons that we should further consider 
the decomposed P-wave channel. One is the decompo- 
sition in Eq.([T]) is not exact, although it holds within a 
few percent, hence we need to study whether this small 
error will be enlarged when comparing with experimen- 
tal data. The other reason is that rg and ri vary with 
different center-of-mass energies or different experimen- 
tal cuts introduced in experiments, thus the two LDMEs 
Mq_^^ and M/^^^ cannot be universally used. Regarding 
this point, we find the changes of rg and ri are not large 
in different cases (see Table |T|. As a result, M^^^ and 
M]'^^^ extracted from the CDF data can be approximately 
used to predict other experimental results. But this can 
also cause some errors. A convenient method to cover all 
these theoretical uncertainties is fitting the experimental 
data using three independent LDMEs. As pointed out 
above, data with pT < 7 GeV may not be well explained 
by the fixed order perturbative QCD calculations, so in 
the fit we still choose p™' = 7 GeV, which is safer for the 
application of perturbative QCD. 

For the J/i/", by minimizing we get 

01 = {O'^/'f'iV^^)) = 15.7 X 10-2 GeV3(±129%), 

02 = (O'^/^Csf^)) = -1.18 X 10-2 GeV3(±249%), (8) 



O. = 



= -2.28 X 10-2 GeV^(±239%). 



These three LDMEs are unphysically determined, which 
is reflected by the large relative errors shown in the end 
of each expressions. Nevertheless, it does not matter 
because we can find some linear combinations of them, 
which are physically determined and have small uncer- 
tainties. Define the correlation matrix C 



1 d^x' 



2 dadO, 



(9) 



at the central value points, we have 



C 



0.041 
-0.0060 
-0.011 



-0.0060 -0.011 
0.00087 0.0016 
0.0016 0.0030 



(10) 



The eigenvalues with corresponding eigenvectors 
of C are then 
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= 4.5 X 10-2 


^1 


= (0.96,-0.14,-0.26) 
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= 1.2 X 10-*^ 




= (0.29,0.31,0.91) 


A3 


= 9.2 X 10-^ 


^3 


= (0.047,0.94,-0.33). 



The LDMEs corresponding to the eigenvectors are 



We have not considered the 5 — - — term in the expansion 



in Eq. [S] which is important in the region of ^ , 
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where we denote matrix 
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(13) 



Inserting Eqs.(|5|) and ([TT|) into Ea.([T^. we have 

Ai = 15.8 X 10^2 Qgv^ (±134%), 
A2 = 2.11 X 10^2 Qgv^ (±5.13%), 



(14) 



0.39 X 10^2 GeV^ (±2.45%). 



It can be seen that A2 and A3 are well constrained in this 
fit, while Ai is badly determined which contains all un- 
physical information in Eq.®. Using A^, the differential 
cross section can be expressed as 



do- = ^ AcTiOi = ^ QiAi 



with = da V-\ 



(15) 
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where d&i denote the corresponding short-distance coef- 
ficients. With its large value and large uncertainty, Ai 
may damage the theoretical results if its coefficient ai is 
not very small. Fortunately, with the CDF condition, we 
find contributions of Ai, , are less than four percent 
for all regions of 7 GeV < pr < 20 GeV. 

In the above treatment, the LDMEs defined in Eq.([2]) 
correspond to vectors \^Mo = (0.25,0,0.97) and \^Afi = 
(0,0.87,-0.48), where we have normalized the vectors. 



We find ^ 



Mo 



^9 and \^ 



Ml 



It means M, 



0,ro 



and M-^ are approximately equivalent to the two well 
constrained ones A2 and A3 respectively. As a result, 
if the badly determined Ai is not important, results of 
using two LDMEs {MqI^ and m(^^) and using three 
LDMEs (Ai, A2 and A3) should be approximately the 
same. 

Comparisons of predictions between using two LDMEs 
and using three LDMEs are shown in FigUand FiglS]for 
the measured CMS [34] and LIICb[35] data respectively. 

For the CMS condition {^/S = 7 TeV and \yj/^\ < 
2.4), we find from Fig|3]that the two methods give almost 
indistinguishable central values and error bars. This is 
understood as ro,i for CMS only have small differences 
from that for CDF, where the LDMEs are extracted. In 
this case, ai in Eq. (|15p is much smaller than 02 and 
as, therefore the contribution of Ai is ignorable although 
it has large uncertainty. We find that the theoretical 
predictions are in good agreement with the CMS data in 
a very wide range of pT ■ 

For the LHCb condition {VS = 7 TeV and 2.5 < 
J/j/i/i < 4), we find from FiglS] that although the two 
methods give the same central values, the method using 
three LDMEs have larger errors when pT > 9 GeV. The 
reason is the influence of a relatively large difference of 
ri between LHCb and CDF (about 18%) on the uncer- 
tainty in the method of using three LDMEs is enhanced 
by the large error of Ai . On the other hand, the rela- 
tively large difference of ri may give a chance to extract 



FIG. 4: (Color online.) Transverse momentum distributions 
of prompt J/tp production at the LHC compared with the 
CMS data for \/^ = 7 TeV and \yj/^\ < 2.4. The CMS data 
are taken from Ref.[33. The two methods give almost the 
same predictions. 
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FIG. 5: (Color online.) Transverse momentum distributions 
of prompt J/^/> production at the LHC compared with the 
LHCb data for = 7 TeV and 2.5 < t/j/^ < 4. The LHCb 
data are taken from Ref . [35l| . The two methods give the same 
predictions for the central values. 



all three LDMEs with small uncertainties when experi- 
mental data at LHCb are adequate enough. Anyway, it 
can be seen from Fig|5]that our predictions give a good 
description for the LHCb data. 

In short, the methods of using two LDMEs and us- 
ing three LDMEs are consistent in giving predictions 
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in the present situation, when only two independent 
LDMEs can be weU constrained. The method of using 
two LDMEs have advantages of simple formalism and in- 
tuitive physical implication, as they approximately rep- 
resent the (for MqI^) and p^"^ (for m(^^) behaviors 
of the cross section, but it needs to consider uncertain- 
ties originated from the decomposition Eq.([T|) and the 
differences of ro.i additionally. On the other hand, the 
method of using three LDMEs can systematically treat 
all uncertainties but with a more complicated form, with 
which it may not be easy to see the physical meaning 
directly. 

Within the method of using two LDMEs, whether a 
good prediction can be achieved is under control from 
the differences of ro.i between conditions under which we 
make predictions and conditions on which the LDMEs 
are extracted. Because the decomposition in Eq.dl]) is 
good in the cases discussed in this work (see Table H]), we 
expect there is no large uncertainty from it. 



IV. PREDICTIONS FOR LHC AND RHIC 

We compare our predictions of J /ip prompt produc- 
tion at the LHC with new LHC data in FiglHl The data 
of ALICE, ATLAS and LHCb Collaborations are taken 
from a recent meeting at CERN [s^ jwhile data of CMS 
Collaboration are taken from Ref. [33 . Besides statis- 
tical and systematic errors, comparable variations from 
spin- alignment uncertainty are also considered in data of 
ALICE, ATLAS and CMS Cohaborations. Errors from 
spin-alignment are dominant for most pt points, therefor, 
more studies on polarizations are needed in the future. 
On the theoretical side, we use the method of using two 
LDMEs as discussed in previous sections. It can be found 
that our predictions are in good agreement with all data 
on the whole. Specifically, from the comparison with 
the LHCb data, we find predicted cross sections become 
declining relative to data as yj/^ becomes larger. This 
phenomenon, however, can be understood easily because 
ri tends to be far away from —0.56 when yj/^ becomes 
larger (see Table H]). On the other hand, as mentioned 
in the last section, the relative large difference of ri may 
give a chance to extract all three LDMEs when LHCb 
has enough data. 

Data at large pt are very important because they may 
distinguish between different models. Recently, both 
ATLAS Hi] and CMSfli] Collaborations have released 
their data of prompt J/ip production for pT as large 
as 70 GeV. Comparisons with our predictions (with 
the same input parameters as in Refs.[ll. [l7j) are shown 
in FiglTl where it is found that all data are located 
within predicted uncertainty bound (a factor of two). 
We fit the CO LDMEs using the Tevatron data with 
7 GeV < Pt < 20 GeV and give a very good prediction 
for the LHC data up to pt = 70 GeV. This is a nontriv- 
ial test for the universality of CO LDMEs. Note that it 
is certainly needed to extract the CO LDMEs from these 



large px data when data are adequate enough. 

Our predictions for tp' prompt production at the LHC 
compared with CMS datafs^ and LHCb dataji^l are 
shown in FiglH] The predictions are in good agreement 
with CMS data. For the LHCb, because the data in- 
clude also B decay contributions, we can not compare 
with them directly, but a consistence between data and 
prediction can still be found. 

We also give predictions for J/4> and -0' productions 
at RHIC in FiglHl It is found that the predictions are in 
good agreement with the data. 



V. COMPARISON WITH RELATED WORK 

Soon after this work was presented in a meeting[l|, 
another talk[42|] (see also [i^) appeared, in which a full 
NLO QCD correction to direct J/ip production was also 
reported. They did not consider feeddown contribu- 
tions of '0(25') and Xcj to J/t/j production, but jointly 
fit the Tevatron data and HERA data for J/ip produc- 
tion (Tevatron data with p^"' = 3 GeV and HERA data 
with p™* — 1 GeV). It is encouraging that, for all 
short-distance coefficients in J/tp direct production at 
the Tevatron, results in our two groups consistent with 
each other. 

However, the results of extracted LDMEs are signifi- 
cantly different. Specifically, they get j42] 

(O'^/'^CSq^)) = (4.76 ± 0.71) X 10-2 GeV^ 
(©■^/'/'(^'fl)) = (0.265 ± 0.091) X 10-2 GeV^ (16) 
(C)^/'/'(3p^8l)) = (-1.32 ± 0.35) X 10-2 GcV^ 

Inserting them into Eq.((T]), we get 

Mg-^^lf = 2.47 X 10-2 GeV^ 

M^-^^f = 0.594 X 10-2 GeV^ (17) 

which are much different from our results in Table HIl The 
authors of Ref. [42] also pointed out that Mq^^ and M^^^ 
are not precisely corresponding to the well constrained 
eigenvectors ^2 and in Eq. pTj) . but also mixed with 
^^1, thus in our fit there are very large uncertainties in 
LDMEs. 

First of all, we note that a small mixing with is 
not so terrible. If we can expect that the physical LDME 
corresponding to is not much larger than that cor- 
responding to ^2 and \^3, then the error caused by the 
mixing is just as large as the size of mixing, a few per- 
cents in our case. When the decomposition of Eq. [1] 
holds very well, there will be a LDME which can only 
be badly constrained. The fitted value of a badly con- 
strained LDME is always much larger than its real value 
because of stochastic effect, which explains the fact that 
LDME corresponding to is much larger than that 
corresponding to ^2 and ^3 in Eq. [TH 
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FIG. 6: (Color online.) Transverse momentum distributions of prompt J/^p production at the LHC compared with the new 
data of ALICE, ATLAS, CMS and LHCb Collaborations for ^ = 7 TeV. The LHC data are taken from Ref.lli,!!^. 
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To clarify the discrepancy between Eg [Tfl and Table 
ini we do a similar fit as authors in Ref. |42[ did: using 
three LDMEs to fit the Tevatron data with pf^ = 7 GeV 
without considering feeddown contributions. Wc then get 

Mq^/^ = 8.54 X 10^2 GeV^ (±12%), 

Mi^/f = 0.167 X 10^2 GeV^ (±63%). (18) 

Comparing this result with that using two LDMEs to do 
the fit without considering feeddown contributions 

Mq^/^ = 8.92 X 10"^ GeV^ (±4.4%), 

m//^ = 0.126 X IQ-^ GeV^ (±18%), (19) 

we find the two methods give very similar Mq(,^ and 

^lit- Comparing Eq.dlll) with Table HI we find the 

feeddown contributions change M^l^ a little but reduce 

Mi^/f by a factor of 2. 

We conclude that, even without subtracting feeddown 
contributions, results of only fitting Tevatron data with 
p™* = 7 GeV in Ea. ([T5)) are still significantly different 

from that in Eq. pT)) . Specifically, Mq(,^ is well con- 
strained in both Eq. pTl) and Ea.([T5|. but the central 
value is much different. The difference, as short-distance 
coefficients are the same and the same fit method is used, 
must be ascribed to different treatments for experimental 
data in the fits. In our opinion, data for pT > 3 GeV 
at the Tevatron and pT > 1 GeV at HERA can not 
be described consistently by the fixed order perturbative 
NRQCD. The inconsistence may imply that the fixed or- 
der perturbative calculation can not describe the data in 
small pt region (3 GeV < pr < 7 GeV for Tevatron and 
PT ^ 1 GeV for HERA). 

Besides, it will be interesting to see if the result given 
in Refs.[43.l43| can describe the large J/ip production 



cross sections (say 20 GeV < pt < 70 GeV) observed 
very recently at the LHC, since the large pt data provide 
a very important test for the LDMEs. 



VI. SUMMARY 

In summary, in this work we calculate the J/ ip and ip' 
prompt production at the Tevatron, LHC, and RHIC at 
0{a'^v'^), including all CS, CO, and feeddown contribu- 
tions. A large K factor of P-wave CO channels at high 
PT results in two linearly combined LDMEs Mq^^ and 

Afi^/f , which can be extracted at NLO from the Teva- 
tron data. We argue that NLO result is necessary and 
essential to give a good description for J/ip production, 
because the NNLO CS contributions are unlikely to be so 
important as to substantially enhance the cross sections 
at large pt- We also compare the method of using two 
LDMEs with that using three LDMEs, and find these two 
methods can give consistent predictions in the present 
situation. For rg^i, which appear in two combinations 
of LDMEs and are related to the short-distance coeffi- 
cients depending on given experimental conditions (e.g., 
the beam energy, the rapidity values,...), when the differ- 
ences of ro,i between the experiment in which the LDMEs 
are extracted, and the experiment for which the predic- 
tion is made, are small, the two methods give almost the 
same predictions with only small errors. Whereas when 
the differences are large, predictions of both of the two 
methods will have large uncertainties. Our theoretical 
predictions are in good agreement with the newly mea- 
sured LHC data and RHIC data for both J/ijj and ip' 
prompt production, which implies that the universality 
of CO LDMEs may hold approximately in charmonium 
hadroproduction. However, if one uses fixed order per- 
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FIG. 8: (Color online.) Transverse momentum distributions of prompt tp' production at the LHC. The CMS data are taken 
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turbative calculation to describe data in the small px 
region, we find the universality of color-octet matrix ele- 
ments may be broken. Our work provides a new test for 
the universality of color-octet matrix elements, and the 
color-octet mechanism in general. 
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